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Abstract 



h^- When a linear system Ax = y is solved by means of iterative meth- 

^~[ ods (mainly CG and GMRES) and the convergence rate is slow, one 

<"^ may consider a preconditioner P and move to the preconditioned system 

"ti P^ 1 Ax = P _1 j/. The use of such preconditioner changes the spectrum of 

a the matrix defining the system and could result into a great acceleration 
of the convergence rate. The construction of optimal rank precondition- 
ers is strongly related to the possibility of splitting AasA = P + R + E, 

y—( where E is a small perturbation and R is of low rank [Tyr96]. In the 

^ present work we extend the black-dot algorithm for the computation of 

C^) such splitting for P circulant (see [OT06]), to the case where P is in j^, 

^O for several known low-complexity matrix algebras j4 ■ The algorithm so 

'-fi obtained is particularly efficient when A is Toeplitz plus Hankel like. We 

^— ' finally discuss in detail the existence and the properties of the decomposi- 

■^- tion A = P + R + E when A is Toeplitz, also extending to the <y3-circulant 

f"*) and Hartley- type cases some results previously known for P circulant. 

m 



1 Introduction 

In this work we consider a new approach to the construction of preconditioning 
matrices P for the solution of linear systems Ax = y. We call this kind of pre- 
conditioners optimal rank because they are produced trying to force the rank 
of A — P to be as small as possible. Optimal rank circulant preconditioners P 
were initially proposed for Toeplitz systems by Tyrtyshnikov et al. in [OT06], 
[ZOT06]. Here we basically extend them to several known low complexity ma- 
trix algebras [DZ01], [BD03], [DFZ06], [Di 09] and to Toeplitz plus Hankel like 
matrices. 

1.1 Notations 

We use some standard notations, which are briefly described here. 



With M„ we denote the Hilbert space ofnxn matrices whose entries are 
in the complex field C, with U n the group of unitary n x n matrices and with 
M+ the cone of positive semi-definite n x n matrices. We use also the notation 
A > for an element of M+ and A > for those matrices which are strictly 
positive definite. Given A e M„, a (A) and \t(A) are the spectrum and the i-th 
eigenvalue of A, respectively. 

The square bracket [ • , • ] : M n x M n — > M„ denotes the commutator 

4,B4 [A, B] = AB- BA 

also denoted by \7 a(B) = [B,A]. The round bracket ( • , • ) denotes the standard 
scalar product on C" (sometimes it can denote the scalar product on different 
Hilbert spaces; it will be clear from the context). 

The symbol e^ is used for the i-th canonical vector (e;)/j = 1 if k = % and 
(ei)fe = otherwise. 

Finally given a matrix W £ M n we shall use the symbol j4{W) for the 
algebra generated by W, namely the closed set 

j4{W) — {p(W) | p polynomials} 

2 Low complexity matrix algebras 

Suppose we are given a unitary matrix U € U n . Then it can be naturally defined 
the algebra of normal matrices 

j* = sdU = {U<&a,g(Oi,...,d n )U* | Bi e C} 

also called algebra of matrices simultaneously diagonalized by a unitary trans- 
form or briefly sd U algebra. For an element A £ j4 , the complexity of the prod- 
ucts A x vector and A^ 1 x vector depends only on the complexity of U x vector 
and U* x vector. Therefore we say that a sd U algebra j4 is of low complex- 
ity if both these products are computable with less than 0{n 2 ) operations, in 
particular if they can be performed with O(nlogn) floats. 

Notice that, whenever W € M„ is diagonalized by M £ M„, for any A <G 
j4iyV) we have A = M diag(Ai(A), . . . , A„(A))M _1 . Moreover if (Ax) is a 
family of mutually commuting matrices, it is known that (Ax) admits a common 
Schur basis. These two facts together imply that the algebra jrf(N) generated 
by a normal n x n matrix JV, must satisfy j4(N) C sd U for a unitary matrix 
[/£ U„. The inclusion can be proper, and precisely it is an identity if and only 
if N is non-derogatory 1 [DZ01]. Therefore the non-derogatorycity hypothesis 
on the matrix W leads to the following further characterization 



jrf(W) =kcrV 



w 



1 A matrix A S M n is said to be non-derogatory if dcg(p) > n, for any polynomial p such 
that p(A) = 0, or, equivalently, if the geometric multiplicity of any eigenvalue of A is one. 
We make often use of matrices that are both normal and non-derogatory, which therefore are 
those normal matrices that have pairwise different eigenvalues. 



and if W is normal we also have j#(W) = sd U, for an U £ U n . 

The best known low complexity sd U algebras are commonly divided into 
three classes: </?-circulants, Trigonometric and Hartley-type. Some specific 
choices in such classes have been used successfully to solve linear algebra prob- 
lems involving Toeplitz matrices or, more generally, structured matrices related 
to shift invariance of the mathematical model considered (see f.i. [CN96] and 
references therein) . Anyway, depending on the problem, any algebra among the 
three families could find a potential application as a preconditioner. 



2.1 (/9-circulant algebras 

Let us consider the matrix 

/ 1 



n ip = 



1 



w 



if£ 



(1) 



/ 



which is the basis for the definition of the family of <p-circulant algebras 

Definition 2.1 Given <p £ C, the algebra G v = j4{U v ) generated by 1J V is 
called ip-circulant algebra. 

Note that ejll^ = e~[ +l for k = 0, 1, . . . , n — 1, which implies that G v is a 
1-space [DZ01] or, in other words, that any matrix C £ G v is uniquely defined 
by its first row. It is now natural to introduce the operator 



C„ 



-V" 



X I 7 0(n I X I 



which maps x £ C n into the matrix C v (x) £ C v whose first row is x T . 
It is not difficult to observe that the matrix 



F, 



1 



f 



ip^uj 11 



(2) 



jn V" / t,j=0,l,...,n-l 

diagonalizes the algebra C v , namely that 

G v = {F v diag^i, . . . , O^F- 1 \6iGC}. 

Moreover F v £ U n if and only if \ip\ = 1 (cf. [Dav79]). From now on we assume 
that if, defining G v , has modulus one, unless otherwise specified. Nevertheless 
we underline that several formulas that we obtain could be adapted to the case 
of a generic complex if. The choice if = 1 gives rise to the well known circulant 
algebra 6, diagonalized by the Fourier matrix 

(3) 



F = 

' k,h=0,l,...,n-l 

which is naturally related to F v up to the diagonal scaling 



F v =diag(l,<£' 



>V 



)F = A V F. 



We shall use simply the symbol II for 77i. Another very popular choice is 
(p = — 1 which defines the so called skew-circulant algebra C_i. Both C and C_i 
and their applications has been widely studied, see for instance [Huc92]. Since 
U v is normal and non-derogatory we have G v = ker Vr/ . 

2.2 Trigonometric algebras 

There are sixteen different trigonometric algebras presently known [Di 09] , [BD95] . 
Eight of them are diagonalized by a discrete sine- type transform, the other eight 
by a cosine-type one. We start again by introducing a family of matrices 



■"■/i ^MMl,M2,^3,A»4) 



( Ml 


/'2 






\ 






1 





1 






, M = 


/ Ml \ 

M2 
M3 






1 


u 


1 




V M4 / 


V 






M3 


M4 ) 







(4) 



and defining T M = ^{X^). Sixteen different choices for the vector (i£l 4 give 
rise to the sixteen trigonometric algebras (see f.i. appendix 1 in [CSST08]). We 
list in Table 1 such particular values for fi naming the corresponding trigono- 
metric algebras T M also DST or DCT so to recall that they are diagonalized 
by a discrete sine transform or a discrete cosine transform, respectively. It is 

Table 1 : Sixteen choices for ^t S K 4 and respective trigonometric algebras. 







M3 = 2 
M4 = 


M3 = 1 
M4 = 


M3 = 1 
M4 = 1 


M3 = 1 
M4 = -1 


Mi = 


= 0, M 2 = 2 


DCT X 


DCT 3 


DCT 5 


DCT 7 


Mi = 


= 0, M2 = 1 


DST 3 


DST X 


DST 7 


DST b 


Mi = 


= 1, M2 = 1 


DCT 6 


DCT 8 


DCT 2 


DCT 4 


Mi = 


-1, M2 = 1 


DST 8 


DSTq 


DST A 


DST 2 



not difficult to observe that for such choices of /i the matrix X^ is normal and 
non-derogatory, since /12 ^ 0, so T M = ker Vx • 

The well known tau-algebra is DSTi = 7^i 10 \ — 1 whose generating ma- 
trix will be denoted simply by X. Such algebra was considered in [BC83] where 
it is defined as the set ofnxn matrices satisfying the cross-sum rule with null 
boundary conditions 



1={AeM n \a tj e 



-i,j 



H+l,j 



H,j-1 



M,j+l 



, i,j = l,...,n 



a-n+ij — ao,j — a«,n+i — flj,o — 0}. 



This is a computationally useful definition and it clearly is nothing but the 
scalar form of our previous characterization of 0" as the kernel of Vx- Other 
trigonometric algebras which attained particular attention are DCT 2 , DST 2 , 
DST3, DST-j mainly because of their applications to image processing [NCT99] 
and displacement decomposition [BD95], [DZ95]. 



As for the 93-circulants, also the trigonometric algebras are 1-spaces [DZ01], 
namely each matrix in T^ is uniquely defined up to its first row. In [DZ95] such 
property is explicitly shown for the more generic set of Hessenberg algebras 
which contains both <p-circulant and trigonometric algebras. As a special case 
of that result we derive the following characterization for T M : consider n — 1 
polynomials $1, . . . , $„_i, each $& defined as the characteristic polynomial of 
the principal submatrix of X^ of order k. Setting $0 = 1 an d 



4 fe) = V2 1 ® 



fc-1 



[X, 



/,: = !,. 



then T„ = Span(A^ , 



, X, 



' n ) and the 1-space property eJX^ 



So, as for the circulant case, we introduce in a natural way the operator 



el follows. 



T„ 



X !->■ Tp(x) 



which maps x E C n into the matrix t^{x) £ 7^ whose first row is x T . 

Since X^ is non-derogatory and normal as well, there exists a matrix U^ € U„ 
such that T^ = sd l/ M . Using the symbol S^ (C„) for the matrix U^ diagonalizing 
T p when the choice of /j, gives rise to a DST (DCT), we list S^, C M and the 
eigenvalues of X^ in Tables 2,3. We finally underline that all of them satisfy 
the low complexity property (see [SGP + 95], [SPS + 96] and references therein). 

2.3 Hartley-type algebras 

In [BD03] eight different unitary matrices Hi are introduced and eight Hartley- 
type algebras !Hj = sd Hi are defined as the set of matrices simultaneously 



Table 2: Discrete sine transform S u and the eigenvalues Ai. of X u , k, h = 0, 



,n — 1 





(S^kh 


-^fe(^) 


DST X 


Hn(7- 1 \\(h 1 \\ 7F 


(fc+lW 


n + 1 


n+1 


DST 2 


sin(fc + !)(/>+ ^ 


(fc+lW 

cos 

n 


DST 3 


sin(fc + i)(^ + l)^ 


cos — 

n 


DST 4 


sin( fc+ i)(/,+ i)^ 


cos — 

n 


DST 5 


7T 
"inf7" 1 1 )(h 1 1 "l 


(k + l)7T 


n +7. 


n+i 


DST 6 


HnO I 1 W 7) I I 7F 


(fc+lW 


""K'x ' 1 M rt ' ) , 1 

V 2 / n+ 5 


~" n+i 


DST 7 


"in 1 /- 1 \ (h \ A\ 


(&+£)tt 


V 2 7 n+i 


n+i 


DST S 




(&+£)* 


oin^ 1 2 j^i 1 2 J n _i 


n- 2 



Table 3: Discrete cosine transform C M and the eigenvalues Afe of X„, k, h = 0, 



,n- 1 





{C^)kh 




Afc(-X^) 


DCT X 


cos kh - 

n — 1 


kit 


n — 1 


DCT 2 


COS K n, H — — 
V 2,/n 


fc7T 

cos — 
n 


DCT 3 


COS K + - 1 ft — 


(fc+^TT 

cos — 

n 


DCTi 


cos (fc +!)(.+ !) 


7T 

n 


(fc+i)7r 

cos — 

n 


DCT 5 


cos fcft t- 


kit 
cos t- 


DCT 6 


, /.. 1\ ^ 




kix 

cos =- 

n- 2 


V 2 / n- i 


DCT 7 


cos \ k H — n, =- 

V 2 7 n-^ 


COS =j — 

n- 2 


DCT 8 




7T 


cos -. — 


cos^+ 2 j^+ 2 j w 


+ 1 



diagonalized by such Hi. The Hartley algebra "K\ was introduced in [BF93] and 
the well known Hartley matrix H\, which diagonalizes it, is defined as follows 



Hi 



1 



\ n 






jj— 0,l,...,n— 1 



The multiplication _ffi x vector can be performed with 0(n log n) operations, 
nonetheless in [BD03] is shown that the same low complexity property holds for 
all the Hi, i = 1, . . . , 8. Let us introduce in detail the Hartley- type transforms 
and the corresponding algebras. Consider the matrices 



K = H 2 = 



G = H* = 



1 



cas 



cas 



2ttz(2j + 1) 



27r(2i + l)(2j + l) 
2n 



ij.— 0,1,. ...n— 1 



ij,— 0,l,...,n— 1 



where casx = cos a; + sinir. Both of them are orthonormal matrices. Further- 
more consider the two sparse n x n matrices 



/ V2 



Ei = 



V2 



\ 

I J 

-J I J 



E 2 = 




where the presence of the central row and column depends on the oddness of n. 
Consider finally the four subspaces of ±l-circulant algebras 



j ±i 



{Cee ±1 \C T = C}, 



psk 



±1 



{Ce e±i \C T =-C}. 



the first two being a proper subalgebra of 6±i. 

The eight Hartley-type transforms and algebras are defined by the identities 
in Table 4. Note that only the algebra § — sdG is not a 1-space, in fact its 
elements are not defined uniquely by the first row (for more details see [DZ01], 
[BD03]). 

Table 4: Definitions for the eight Hartley-type algebras M; 



5Ci = 'K = sd H = sd H x 



jne s 



5i 2 = % = sdK = sdH 2 = Cij + Ji7_ iCi^ 



3i 3 = g = sdG = sdH 3 = e s : 



je s : 



IK4 


= %■■ 


= sdK T = sd H 4 = G S + je sk 


3^-5 


= V = 


= sd(K T E 1 ) = sdH 5 = e s + je s 


!K 6 


= A* = 


= sd(GE 2 ) = sdH 6 = eij + jeii 


3^7 


= a = 


= sd{HE[) = sdH 7 = e s + jne s 



3-Cs 



f3 = sd{KEj) = sdiJ 8 = Cij + JiJ-iCix 



Unlike 93-circulant and trigonometric algebras, it is not so clear that Hartley- 
type algebras can be introduced as the algebras generated by matrices whose 
structure is predictable for all n. However, since all of them are algebras of 
normal matrices simultaneously diagonalized by a unitary transform, there exist 
non-derogatory matrices W{ such that 3-Q = j4{W%) = kerV^v Nevertheless 
let us note that the C^ part of Jf^ is the subalgebra j4(Y±\) generated by the 
following derogatory matrix 



Y : , 



n a + n' 



( 

1 

W 



1 



Noting that 



f^h p f v = o v = p 1/n n = ^/"diag^ | % = o,. 



*\ 



1 
/ 



,n-l) 



(5) 



and that U±\ + ZL^ = II ±\ + iTJj = 2iRi7±i, 2 we easily obtain an explicit 
formula for the eigenvalues of Y\ and Y_ 1 , namely 



for k = 0, 



Afc(Fi) = 
. . , n — 1. 



2 cos 



(£) 



A fe (F_i 



2 cos 



(2fc+l)7T 



(6) 



2 The real part (Hermitian part) of a matrix X is the Hermitian matrix 5RX = ■= (X + X*) 



3 Optimal rank preconditioning 

Let us consider an n x n linear system 

Ax = y, AeM„, xjeC™ (7) 

which should be solved by some Krylov subspace iterative methods (CG and 
GMRES are good examples). When the convergence rate of such methods 
is low, one may consider a suitable preconditioning matrix P, switch to the 
preconditioned linear system 

p- l Ax = p- 1 y, (8) 

and apply to this new system the iterative methods. Clearly except for trivial 
cases, the spectrum of the matrix defining the system (8) is different from the 
original one, while the solution x maintains unchanged. The introduction of the 
preconditioner P could lead to a substantial improvement of the convergence 
rate, provided that P satisfies some "good" properties, that we summarize into 
the following two 

1. The condition number k(P~ 1 A) is uniformly bounded in n 

2. The spectrum of P~ l A has a cluster around 1. 

When A > (the general case needs some further hypothesis cf. [TZY97]), the 
first property leads to the linear convergence of the methods, independently on 
the dimension n of the problem. The second one is related with the super-linear 
convergence of the methods and, fixed an e > 0, it is the same as requiring that 
the following splitting for the matrix A in (7) 

A= P+R+E 

holds, being E a small perturbation, \\E\\ < e, and R a low-rank matrix, namely 
rank R is o(n) [Tyr96] , [OT06] . It is clear that the cluster oia{P~ l A) is a proper 
cluster whenever rank R = r e (n) is uniformly bounded with respect n, in fact 
r s (n) is exactly the number of eigenvalues of P _1 A which are outside a ball 
of radius e around 1. Moreover we can heuristically affirm that "the smaller is 
rankR, the smaller is the cluster of a(P~ 1 A) and the better is the preconditioner 
P". 

In the following we consider low complexity algebras of matrices simultane- 
ously diagonalized by a unitary fast transform U, or rather closed sets of the 
form 

j* = sdU = {Udi&g(0i,...,6 n )U* | 6i€ C} 

where U x vector and U* x vector require O(nlogn) operations. 

Our considerations about property 2 above suggest the definition of an op- 
timal rank matrix algebra preconditioner for a given linear system (7) . 



Definition 3.1 Given an invertible n x n matrix A and a matrix algebra j4 , 
we call optimal rank preconditioner in j4 for A any matrix 

A = argmin{rank(^ -P-E)\PGjs/, \\E\\ < s} . 

In [OT06] it is stressed that in order to construct such preconditioner for A 
one should solve the following 

Problem 1 Given A e M„, given an sdU algebra jtf and given e > 0, find 
A = P + R, with P e j4 and ||yl — *4|| < e such that ranki? is as small as 
possible. 

Observe that for any given M € sd U there exists a diagonal matrix D such 
that M = UDU*. Therefore if the norm considered in Problem 1 || • || is 
unitarily invariant 

||A-.A|| = \\U*AU -D-R\\, rank R = rank R . 

In other words we can split problem 1 into the following problems 2 and 3 and 
calculate the minimum over the algebra of diagonal n x n matrices @ 

Problem 2 Given A e M„ and e > 0, find A = D + R, with D e © and 
111^4 — A\ || < e such that ranki? is minimum, for a chosen unitarily invariant 
norm |j| • |||. 

Problem 3 Given A <G M n and given an sd U algebra j4 = sd U, compute the 
image U*AU. 

If jrf = sdU is of low complexity, the computation of U* AU requires an 
amount of 0(n 2 logn) operations which is not acceptable. For this reason we 
have posed Problem 3. However, we shall see that in solving Problem 3 it is not 
necessary to compute all the entries of the image matrix U*AU, it is instead 
enough to have an algorithm that computes any prescribed entry (U*AU)ij in 
a fast way, by a number of operations independent of the matrix size. 

Problem 2 can be approached as in the circulant-Toeplitz case [OT06], i.e. 
by means of the black-dot algorithm, an ad-hoc version of the incomplete cross 
algorithm [GTZ97], [TyrOO], [BebOO]. Given an algebra j4 = sdU = j4(W) (W 
non-derogatory), in order to apply the black-dot algorithm we need to know 
some elements (U*AU)ij, where the pair (i,j) belongs to a certain set of in- 
dices Q. Typical choices for Q give rise to a method whose complexity can be 
estimated with 0{nr e (n) 2 ) [OT06]. When the coefficients matrix is positive 
definite we expect that the preconditioner A computed by such algorithm is 
positive definite as well. Actually such property is not always ensured, in fact 
it may happen that some entry of the computed diagonal matrix D is negative. 
However, in typical cases - see [OT06] and the proposition here below - it is 
possible to ensure the positive definitess of A by applying a low-rank correction 
to the computed D, i.e. modifying a small number of its diagonal entries. 



Proposition 3.2 Let A = P + R + E with ranki? = r < n and \\E\\ < e. If 
A > el and R > then at least n — r eigenvalues of P are positive. 



Proof. Weyl's inequalities for the Hermitian eigenvalues problem Z = X + Y 
give us the following 

K+j-i(Z) < \i{X) + Xj{Y) , l<i + j-l<n, 
\i(X) + \ n {Y) < \i{Z) <\ i {X) + X 1 {Y), 1 < i < n , 

where the eigenvalues of a Hermitian matrix are supposed in decreasing order 
(i.e. Xi > Ai + i). Since r < n, we have A r +fc(i?) = for k = 1, . . . ,n — r. Thus 
Weyl's inequalities applied to A = P + R + E give us 

X t+r (A) = max X i+k (A) < X t (P + E) < Xi(A), l<i<n-k. 

k— r,....n — 1 

Therefore we know that all the first n — r eigenvalues of P + E are greater or 
equal to e, thus n — r eigenvalues of P are positive, since A > el. □ 

In the next section we discuss Problem 3 and propose a method for com- 
puting any element (U*AU)ij at a very low cost, after a preprocessing phase 
of complexity O(nlogn). Such method works well if the displacement rank of 
A, the matrix defining our system (7), with respect to the chosen algebra s4 ', is 
sufficiently small. More precisely let us introduce the following 

Definition 3.3 Let j4 = j4{W) be a matrix algebra. We say that a matrix 
A £ M„ almost-belongs to j4 , in symbols A G j4 , if rank([A, W]) is uniformly 
bounded in n. 

Note that by the fundamental theorem of homomorphism we have a canon- 
ical isomorphism Vvi/ between M n/j^ and range Vyi/ which implies that the pre- 
image of A e range Vw is given by the closed set {V~^{A) + A}aejJ- 

The method proposed in the next section for computing cheaply the elements 
(U*AU)ij, works if A e s4 and [A, W] is explicitly known. We underline since 
now that our method works well for all the low complexity sd U algebras j4 
previously presented and for any Toeplitz or Hankel matrix A. 

4 The computation of U*AU 

Let j4 C M„ be an algebra of normal matrices j4 = sdU, U € U(n). Let us 
consider an element W € j4 which is non-derogatory. Clearly j4 = j#(W) = 
kerVw and U*WU £ @. Set D = U*WU. Given a matrix A £M n , we have 

(U*AU)D - D(U*AU) = U* [A, W]U 

thus the off-diagonal elements of U*AU satisfy the identity 

(U AU) %J - Xj{w) _ K{w y *¥>3- W 

10 



Such equalities let us state the theorem below, whose detailed proof can be 
derived also by the observations which follow. 

Theorem 4.1 Let j4 = jtf(W) be a low complexity sd U algebra, Aci and 
rank[A, W] = p (thus p is uniformly bounded in n and Xi(W) — Xk(W) <=> i = 
k). Assume that a(W) is explicitly known. Then, after a preprocessing phase of 
complexity O(nlogn) required for the computation of U*[A,W]U (see below), 
each off-diagonal element (U*AU)ij can be calculated with p+1 multiplications. 

The theorem above gives us a tool for approaching Problem 3, or rather for 
computing the image of a matrix A under a unitary fast transformation. In 
fact, we are requested to compute 2p transformations Xk = U*Xk, y~k = U*yk, 
for k = 1, . . . , p, in a preprocessing phase, where Xk and y^ are the vectors 
defining a dyadic decomposition of [A, W]. Then each off-diagonal entry of 
U*AU is known up to p + 1 multiplicative operations (U*AU)ij = (\j(W) — 
\i(W))~ l "Y%-i{xk)i{yk)j- Finally observe that in order to apply the black-dot 
algorithm we actually do not need the diagonal entries of U* A U [OT06], thus 
Theorem 4.1 is not restrictive in our situation. 

To summarize, in order to compute the optimal rank preconditioner A for A 
into the algebra j4 = sd U = ker Vvi/j we can propose the following DR-scheme, 
where the black-dot algorithm is like a black box which sometime requires an 
entry (U*AU) l3 , i^j. 

DR-scheme 

Assume sdU = kci\7 w and [A, W] — J2k=i x kVk- 

1. (Preprocessing) Compute the 2p fast transforms Xk — U*Xk, yk — U*yk, 
k = l,...,p. 

2. Start the black-dot iterations 

2.1. if (U*AU)ij is required, compute it via the identity 



(TT* ATT\ — ^k=l\ X k)iyyk)j 

{U AUhl ~ X j (W)-X i (W) 
where z denotes the complex conjugate of z. 

3. Proceed with the iterations until convergence, passing through 2.1 if nec- 
essary 

Note that the preprocessing phase requires 0(pn log n) operations, whereas only 
p+1 arithmetic multiplications are needed each time step 2.1 must be performed. 
Despite the somewhat general formulation of the DR-scheme, in the following 
we consider some specific cases, in which A is a Toeplitz or Hankel matrix, 
and discuss how to compute U*AU explicitly, for the low complexity algebras 
presented in Section 2, underlining that for these choices of A, we actually have 
A I j4. 
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4.1 (/9-circulant algebras 

Recall that the generic (p-circulant algebra is defined as 

G v = j#(II v ) = &AF V = kcr Vn v , 

where F v — A V F — diag(l, tpn , . . . , tp~^~)F, \tp\ = 1, and F is the Fourier 
matrix (3). We will make use of the following matrix J, also called reverse 
identity 

( l \ 

1 

J = 

Vi J 

Given p € Z, consider the equivalence class {p mod n} and let [p] n denote its 
unique representative in {0, 1, . . . , n — 1}. 

PROPOSITION 4.2 Let T n = (ti-j)ij be a Toeplitz matrix and U^ the matrix 
(1) generating C v . Then T„ e C v and 

[T n ,II v }=x v e[+e n yT (10) 

being a;^ = {tpti- n -ti,... ,<pt-i - £„_i,0) T and y v = -Jx v . 
Proof. Due to the definition of n v we have the following equality 

n-l 

[T n , n v }ij — 2_^ (U-k{n v )kj — {n v )iktk-j) = ft i _^ +n _^ n — *[i+i]„— j ■ 

fe=0 

Therefore (10) holds for i ^ n — 1 and j ^ 0. In fact for such choices of indexes 
we have i — \j + n — 1]„ = i — (j — 1) = [i + 1]„ — j thus [T„, i7 y ]jj = and 
rank([T„, II V ]) — 2. It is not difficult to observe that (10) also holds for i = n — 1 
and j = 0. □ 

The above Proposition immediately implies that, when A has Toeplitz struc- 
ture, the DR-scheme can be applied to the case U — F v , the unitary matrix 
diagonalizing C v . The Hankel case needs some further observation since the 
rank of [H n , II V ] is not bounded in general and a direct use of the DR-scheme 
would be prohibitive. Given a Hankel matrix H n call T(H n ) the Toeplitz ma- 
trix JH n . Here below we observe that when tp is 1 or — 1, the computation 
of F*H n F v can be brought back to the Toeplitz case F*T(H n )F v , for which, 
instead, the DR-scheme works well. 

Observe that 

/ Tl-2fc-l \ 

A JA V = diag I tp « \k = 0, ...,n-l)J, 



^-l 



therefore when tp e { — 1, 1} we get the equality F*JF V — tp « F* JF. More- 
over, since F — JIIF*, F 2 = J II and F* JF = Q, we also have 



f* jf = (F*nF)(F*y = njn , 
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where Q = diag(l, w, . . . ,u) n 1 ), w = e 2,T1 />\ N ow use the definition of T(H n ) 
to write 

J^ff„i^ - (F;jF„)(F;T(H n )F v ) . 

The formulas obtained so far imply the desired result: 

{F%H n F v )a = ^^V- 4 J"(F;T(iI n )F v ) MW , 
i,j = 0, ...,n- 1, i^j. 

Before proceeding to our discussion for the trigonometric and Hartley cases, 
let us introduce some further useful notation. Given a, b £ C™ with ai = &i, 
we shall denote with T n (a,b) the Toeplitz matrix whose first column is a and 
whose first row is b T . Analogously, given u,v £ C™ such that u n — v\ we shall 
denote with H n (u,v) the Hankel matrix whose first row is u T and whose last 
column is v. By Proposition 4.2 the following formula holds for any a, b £ C" 
with ai = fei 

[T n (o,6),/Z v ] = ^(o,6)- J)K v {a,b) T J, (H) 

where )K v (a, b) is the rank one matrix 



jK v {a 1 b) = (<^J6 — II !p a)e 1 . 



(12) 



4.2 Trigonometric algebras 

As discussed in Section 2 all the sixteen trigonometric algebras T M are generated 
by the matrix X^ (4) , for the sixteen choices of \x £ R 4 shown in Table 1 . Being 
X^ non-derogatory and normal all its eigenvalues are distinct and we know 
them explicitly (Tables 2 and 3). As a consequence, Theorem 4.1 holds for the 
trigonometric algebras and for the set of matrices A £ M„ such that A £ T„. 
Let us show that all the Toeplitz and Hankel matrices belong to such set. 



Given /i £ R , let us split X^ into X^ = X 



Mfj,, where 



X 



( 



\ 



1 



\ 



/ 



M„ 



( Mi 

V 



M2 - 1 



A»3 - 1 M4 / 



For any two vectors a, 6 £ C n with ai = b\ we obviously have [T„(a, b), X^] = 
[T n (a, b),X] + [T n (a, b), M^}. Thus it is possible to prove the result for the tau 
algebra 7 and then for all the other algebras T M . 

PROPOSITION 4.3 Let A £ M n be Toeplitz, Hankel or the sum of them. Then 
A £ T M for any /j, £ C 4 , precisely rank[A, X M ] < 8. 

Proof. Notice that rank(M M ) < 2 for all \x £ C 4 , thus rank(L4, M M ]) < 4. 



Moreover note that given any Toeplitz T = (ti-j)ij or Hankel H = (h 



l i+j )iji 



the 
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matrix T + 77 satisfies the cross-sum rule with non-null boundary conditions, 
that is 

U-i-j + hi-i+j + ti+i-j + hi+i+j — (ti-j+i + hi+j-i + ti-j-i + hi + j + i) = 

for any i,j <G Z. Therefore the boundary conditions are given exactly by the 
border columns and rows of the (n + 2) x (n + 2) matrix embedding the given 
T + 77 nxn matrix and maintaining its same structure. Such columns and rows 
can not be null except for trivial cases. As a consequence, since A G T if and 
only if it satisfies the cross-sum with null boundary conditions, we have that 
rank[T„ + 77„, X] < 4, namely 



[Tn 



H n , X] 



* 




* 




O n -2 




* 




* 



where * are in general non-null entries and O n -2 i s the null matrix of order 
n-2. □ 

Proposition 4.4 Let a, b,c,d G C" with a\ = &i and c„ = d\. Then 

[T n (a,b),X} = 6(a,b)-J6(a,b) T J, (13) 

[H n (c,d),X} = 6{d,Jc)J-J0(d,Jc) T (14) 

where 0{a,b) is the rank two matrix 

8(a 1 b) = e l (IIb) T -{na)e[. (15) 

Proof. Displacement formula (14) for H n (c, d) clearly follows from the previous 
one (13), since H n (c,d) = T n (d, Jc)J and 

[H n (c,d),X] = [T n {d,Jc),X]J = 8{d,Jc)J - JO{d,Jc) T . 

By Proposition 4.3 only the border rows and columns of [T n ,X] are non-null, 
for a Toeplitz matrix T n . Thus we just need to check (13) for such four vectors. 
For instance, noting that X = IIq + ZZg 7 ", we have 

[T n (a,b),X]ei = T„(a,6)e 2 - Xa= 77 T a + b 2 e 1 - Xa = b 2 e x - 17 a 

and analogously 

e{[T n (a,b),X] =b T X-e 2 r T n (a,b) = (II b) T - a 2 ej . 

Therefore (13) holds for the first row and column, if we define 0(a,b) with IIq 
in place of 77. It is not difficult to observe that the same can be said also for 
the the last row and column. Moreover, thanks to arithmetic cancellations for 
the corner positions (l,n) and (n, 1), identity (13) holds for our definition of 
0(a, b), given in terms of 77. □ 
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The above proposition gives us an explicit formula for the displacement rank 
of a Toeplitz or Hankel matrix into the tau-algebra X Then, if S — S^o, 1.1,0) 
is the sine transform which diagonalizes T, we can apply the DR-scheme to the 
case where U = S and A = T n + H n is the sum of any two Toeplitz and Hankel 
matrices. 

Anyway, by Proposition 4.3, the displacement rank of a Toeplitz or Hankel 
matrix into any trigonometric algebra 7^ does not exceed 8. Let us first derive 
an explicit formula for the dyadic decomposition of [T n , X^\ . Since 

Mp = ei(/iiei + (n2 - l)e2) T +e„((M3 - l)e„_i +^ 4 en) T , 
we have 

[T n (o,&),M^] = ci(a,&)f/iiei + (/U 2 - l)e 2 J 

+ c n (a, b) ((fi 3 - l)e„_i + H4e n ) 

- ei f /xici(a, b) + (/i 2 - l)c 2 (a, b)j 

- e„n/X3 - l)c„_i(a, 6) + ^ 4 c„(a, 6) J 

where Cfe(a, 6) is the fc-th column of T n {a, b). The required formula is obtained 
by summing the latter one and (13). A similar computation provides a formula 
for [H n ,Xfj] and thus the DR-scheme can be applied when U = S^,C^ is any 
trigonometric transform and A = T n + H n . 

4.3 Hartley-type algebras 

Let "Kk = S( i Hk denote a generic Hartley-type algebra, we can characterize "Kk 
as the set [Di 00] 

■K k = {A I [A, Y v ] = [A, M fe ] = 0} (16) 

where ip E {1, —1} and Mj, S M„ depend on the Hartley-type algebra. We are 
considering. Despite it is obviously possible to define "Kf, as the set kcrVw for 
some non derogatory matrix W, it is not always easy to find a "simple" W with 
such property. Therefore here we make use of (16) and derive an easy variant 
of the DR-scheme. 

Theorem 4.5 Let X and Y be two distinct normal matrices which commute. 
Let Xi(X), Xi(Y) be the eigenvalues of X, Y corresponding to the same common 
eigenvector. Consider the algebra j4 = {A | [A, X] = [A, Y] = 0}. Then j4 is n 
dimensional if and only if Xi(X) = Xj(X) implies Xi(Y) ^ Xj(Y), for any pair 
of distinct indices (i,j)- 

Proof. Denote with Ai, . . . , A„ and /ii, . . . , \i n the eigenvalues of X and Y, re- 
spectively. Since [X, Y] = we have X — UD\U* and Y — UD^U* , where D\ 
and D^ are the diagonal matrices such that (D\)a = Xi and (D^u = \Xi- 
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First let us prove the implication (<=). Take A e j4 , then U* AUD\ = 
D X U*AU and U*AUD f _ l = D f _ l U*AU. Let (U* AU)^ = %. By writing the first 
relation entry wise we get a,ij(Xj — Aj) = 0, therefore a^ must be zero if A^ 7^ Aj. 
When Xi = Xj we use the second relation obtaining a^dij — jii) = 0, which gives 
us Ojj = due to our hypothesis. As a consequence we have ay = Vi ^ j, 
hence U* AU £ @ and dim.*/ = dimsdC/ = n. 

Vice versa, since any A € sd £/ is an element of j4 and dim j4 = n, we have 
j?/ = sd [/. Now proceed by absurd and assume the claim to be false. Without 
loss of generality suppose Ai = A2 and [i\ = /12, an d consider the matrix 



A = U 



( di b \ 

d 2 



U* = UBU* 



V d n J 



where b ^ 0. We have t/*L4,A]£7 = BD X - D X B = and U*[A,Y]U = 
BD^ — D^B = 0. Therefore A commutes with both X and Y, namely A <G j4 . 
This is impossible since A $ sd U. D 

By the above theorem it is clear how to adapt the DR-scheme for Hartley- 
type algebras 

Assume [A, Y v ] = Y, p s=1 x s y* s and [A, M k ] = YZ=i w * z *s- 

1. (Preprocessing) Compute the 2(p + r) fast transforms x s = H^x s , y s = 
H%y s , w t = fffeW t , 2i = H%z t , s = l,...,p,t = l,...,r. 

2. Start the black-dot iterations 

2.1. if (H£AHk)ij is required, compute it via the identity 



(H* k AH k ) i:j 
if X t (Y v ) ^ X,(Y V ), or 

{ H *k AH k)ij -- 



Aj(^v) ~ Ai(V ¥ ,) 



EI=i(w s )i(i 



otherwise 

3. Proceed with the iterations until convergence, passing through 2.1 if nec- 
essary 

Note that in this modified version the DR-scheme requires one more assumption, 
since the matrix A we are considering should have small displacement rank with 
respect to two different matrices Y v and M k , rather than only one. Of course 
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there are many matrices which satisfy such assumption, however it is not obvious 
if among them there are also Toeplitz or Hankel matrices. The rest of the section 
is devoted to observe that this is true in particular cases. 

First of all, note that the matrices Y v are derogatory and only dim C* of their 
eigenvalues are distinct, in fact Y v is the generator of the family of algebras 

Note that in (6) we derived an explicit formula for the spectrum of Y±\. 

Moreover it is clear from Propositions 4.2 that for any Toeplitz matrix A 
we have rank [A, Y v \ < 4. The same conclusion holds for [A, Y v ] with A Hankel. 
Namely 

PROPOSITION 4.6 Assume \ip\ = 1. Let a, b e C™ such that a\ = b\. Then 

[T n (a, b),Y v ] = Hi v {a, b) - JK v (b, a) T + j(tf< v (b, a) - XC v (a, 6) T ) J 

where jfC v {x, y) is the rank one matrix in (12). 
Let c, d € C™ with c n = d\ . Then 

[H n (c,d),Y v ] =e 1 d T (JII-<pI)-(jn-<pI)del+e n c T (IIJ-tpI)-(nj-tpI)ceT 

Proof. The formula for T n (a,b) immediately follows from (11), noting that 
T„(a,fe) T = T n (b,a) and that [A, B T ] = -[A T ,B] T , V A, B e M„. The second 
one, for H n (c, d), follows from (14) and the decomposition Y v = X + R v , being 
R v = (p(e n e{ + eieX)- In fact H n (c,d)e n = d, H n {c,d)e\ = c, e[H n {c,d) = c T 
and elH n {c,d) = d T . U 

Clearly we can also represent [T n (a, b), Y v ] by means of 0{a, b), in fact 

[T n (a, b),Y v ] = 9(a, b) - J6(a, b) T J + (pj(bej - ei b T ) + (p(ae{ - e ia T )J 

which easily comes from (13). 

Now, concerning M&, let us fix four choices among the eight possible indices 
k, precisely k = 1, 2, 5, 6, and consider the corresponding Hartley-type algebras 
"K = JCi, % = 3^2, rj = J-Cs and /i = "Kq. For such choices we can characterize 
Mfc somehow explicitly, in fact in [DZ01] it is shown that 

M * = '' + (oX>) (17) 

where r(zk) — T (o.i,i,o)( z fc) € M(n — 1) and the vectors Zk & IR™" 1 are, respec- 
tively, 

z\ = |(e 2 -e„_i), Z2 = -|(e2 + e„_i), z 5 = z 6 = . 

As a consequence we get t{z§) — t(zq) = O n -\, and 

t{zi) = \X{I-J), t(z 2 ) = -\X{I + J) 

where X = X(o,i,i,o) € M(n— 1) is defined in (4). Also note that for the algebras 
IK5 = rj and JCg = jjl a more elegant characterization does hold. We state it by 
means of the following 
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Proposition 4.7 Let 77 and /x be the Hartley-type algebras defined in Table 
4. Then 

77 = ker Vy 1+ ./ ,u = ker Vy_ 1+ j 

Proof. Clearly 77 C kerVyj+j, /i C ker Vy_ 1 +j and the equalities hold if and 
only if Yi + J and K_i + J are non-derogatory. If i? 5 is the unitary matrix 
diagonalizing 77, it is not difficult to observe that 



H£JH 5 = 



1 11 



where m = ri/2 if n is even and m = (n + l)/2 otherwise. This remark and 
(6) imply that the eigenvalues Ai(Yi + J) = Ai(Yi) + Aj(J) are all distinct; thus 
Y\ + J is non derogatory, and the thesis follows for 77. In the same way one 
proves the thesis also for /i. □ 

It is now clear that we can apply the modified DR-scheme to a quite general 
class of Toeplitz and Hankel matrices. In fact for any symmetric Toeplitz matrix 
T n and any persymmetric Hankel matrix H n we have [T n ,J] — [H n ,J] = 0, 
therefore [T n ,XJ] = [T n ,X]J and [H n ,XJ] = [H n ,X]J. Using (13) and (14), 
it is now straightforward to derive the formulas for the commutator of T n and H n 
with Mfc, for k = 1, 2, 5, 6, taking into account that the matrix X which appears 
into (17) has order n — 1. This eventually allows us to apply the modified 

DR-scheme to the case A = Toeplitz symmetric + Hankel persymmetric. 

We conclude this section by noting that using only Y v and formula (9) we 
can compute at least 0(n(n — 2)) entries of H*AHi, for any i = 1, . . . , 8. The 
elements • of H* AHi that we can not compute this way, are in the positions 
shown in the figure below (we represent them for n = 5, 6). 



(' 



(' 



\ 



( ' 

V 

/ • 

V . 



. ) 



(v = i) 



(*> = -!) 



5 Algebra-plus-low-rank approximation of a ma- 
trix 

As we underlined the previous section, an optimal rank preconditioner P should 
realize the splitting 

A = P + R + E (18) 
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where ||2?|| < e and ranki? is minimum. When A is Toeplitz, it can be 
shown that optimal or Strang-type preconditioners P, chosen inside suitable 
sd U algebras, realize an analogous decomposition where in general rank R = 
0(e' p )0(n), for a p > [Tyr96], [DZ01], [CN96]. In this section we show that 
for particular classes of matrices A and algebras j4 there exists P e j4 that 
realizes the splitting (18) with ||.E|| < e and ranki? = o(n). 

If T is the unit circle T = {z <G C | \z\ — 1}, let us denote with T ni H n : 
L°°(T, C) —} M„ the Toeplitz and Hankel operators, respectively, which map 
/ e L°°(T, C) into the n x n Toeplitz or Hankel matrices T n (f), H n (f). Finally 
call £(T) the subset of L°°(T,C) of all piecewise holomorphic functions with 
logarithmic singularities, i.e. functions given by an holomorphic function plus 
a function with logarithmic singularities. A generic / G £(T) has the form 



f(z) = g(z) + ^2 ^2 a kh ■ (z - z h ) k log(z - z h ), z e T , 



k=0h=0 



with g holomorphic over a set containing T and Zh € T, h = 0, . . . , q. 

In [ZOT06] it is shown that for an / e £(T) the Toeplitz matrix T„(/) 
admits the decomposition 



T n {f) = P + R + E, 



rank(_R) = O ( log - I log — h log n 



where P is a circulant matrix and E is a small perturbation as usual. The space 
£(T) is a special class of symbol functions that, however, covers all examples 
considered in literature on superlinear preconditioners [OT06]. 

In this section we will show that a splitting analogous to the one in [ZOT06] 
also holds for H n (f) and P chosen inside a generic ^-circulant algebra, (p e T. 
We will explicitly describe such matrix P in several cases, and we will discuss 
also the case of Hartley- type algebras. 

Let us consider a A G C and define the vector 



p(X) = ( 1 A A 2 



A' 1 



It is not difficult to observe that the Toeplitz matrix 

r n (p(A),A 1 -"Jp(A))=p(A)p(A- 1 ) T 

is a rank one matrix. Moreover, by requiring a Toeplitz matrix to be of rank 
one, we observe in fact that 

rankT„(a, b) = 1 «=> 3A | a = p(X) and b = A 1 "" Jp{\) . 

Given a real number A set 

( \ ^ 

A 1 

Z n W =T n (p(\),e 1 ) = 



\ \ n - 



A I ) 
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This is the Toeplitz matrix generated by the symbol 



1 
1 -~Ae ; ' 



&W) = ^-T7iS> A >^ 



By noting that p(X) satisfies 

n v p(\) = Xp(X) + (<p - X n )e n 
we obtain the following identities 

^(p(A),ei) = <pe n ej - n v p(\)el = (A"e„ - Ap(A))e 1 T , 
M v (p(\), X l - n Jp(X)) = ((^A 1 -" - X)p(X) -(if- A")e„) ei T 

which lead to 

Lemma 5.1 For any A <G M there exists P v e C v such that the triangular 
Toeplitz matrix Z n (X) splits into Z n (X) = P v + R, with ranki? = 1. 

Proof. Using (11) and the formulas above we have 

r A™ 

[T n (p(x), ei ),n v ] = \ J -^T n (p(x) 1 x 1 - n jp{x)),n^ 

Notice, in fact, that, after the subtraction between )K V and JJfC J, the term 
<pe n ej is canceled. D 

Observe that we can explicitly write the <p-circulant matrix P v in Lemma 
5.1, in fact it is given by the difference 

A™ 
P v = T n (p(X), ei) - ^—T n (p(X), A 1 -" J P (X)) . (19) 

Hence, using the operator C v : C™ — > C v defined in Section 2.1, we have the 
equality 

Z„(A) - ^(.t v (A)) + R, x v (X) = -^-^JII v p{X) . 

As a consequence, the following characterization holds 
Proposition 5.2 The eigenvalues of C ip (x ip (X)) are 

/27Tfc\ 

(a , k = 0,...,n-l. 

V/" V « / 

Proof. Observe that for a generic f the following formulas hold 

n v p{v) = vp(y) + ((f- v n )e n , Jp{v) = v n - x p{y- x ) . 
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Then note that 



nFlx^X) 



v "' FA V (XJ P (X) + &- A n )ei) = ^^FA^X- 1 ) 



(p — X 



x n 



ip-X r 



l-yA~" 



V 



ip — X 
X 

\ 

I 



-FA v Jn vP (X) 



ip- X n 

( 



1-A^- 1 /" 



V 



\ 

I 



The thesis follows by recalling that, given yeC", the eigenvalues of C v {y) are 
the entries of the vector ^/nFjy. □ 

Consider now the matrix 

K n (X) = Z n {X) + Z n (X) T -I, 
i.e. the well known Kac-Murdock-Szego (KMS) matrix generated by the symbol 

1-A 2 



K\ 



(d) = 2Vi(x(8)-l = J2 xlnleine 



1 - 2A cos 6 + A 2 



(20) 



Observe that when A is real, we have the identity, K n (X) — 2 s RZ n (X) — /. 
It is also easy to show that k\(8) > •*=>■ |A| < 1 or, equivalently, K n (X) > 
«=> |A| < 1. Therefore 

Proposition 5.3 For a given A e R, let if„(A) be the corresponding KMS 
matrix. Then 

1. For any ip <E T there exist Q v € G v and i? of rank 2 such that if„(A) = 

2. Let £ v (\) = (^Jn v + =^) p(X) - 1. Then Q v = C„(^(A)) 

3. The eigenvalues of Q v are /c a (^^), fc = 0, 



, n — 1 



Notice that we also have K n (X) > •<=>■ Q v > 0. 

Proof. By Proposition 5.1 we have the equality K n (X) = (2 s RC v> (x tp (X)) — I) + R 
where R is a rank two Hermitian matrix. Note that the matrix RC v {x V (X)) 
belongs to S„, and this proves (1). Moreover its eigenvalues are the real part of 
the eigenvalues of C v (x v (X)). Therefore from (20) we derive (3). Concerning 
(2) just observe that, for any y <G C", C v (y)* — TpC v (jn v y), and use such 
remark to compute the first row of 2$lC lp (x v (A)) — I. □ 
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It is important to note that the ^-circulant matrix Q v in the previous propo- 
sition is indeed the optimal rank (/5-circulant preconditioner for a KMS matrix. 
This fact can be easily proved by a direct calculation. Just try to impose that 
the difference between K n (X) and a rank one matrix is </?-circulant to reach an 
absurd. Therefore we have an explicit formula for the optimal rank precondi- 
tioner of K n (X). Note that it outperforms, from the clustering point of view, 
any other known preconditioner for a KMS matrix [CN96] , [SE87] . In fact the 
preconditioned matrix Q~ K n [X) has only three distinct eigenvalues. 

Theorem 5.4 Let p, q be two complex valued mutually prime polynomials 
defined on T, such that ^ <z(T), dcgp < degq, and q has all distinct roots. 
Then, for all ip £ T, the lower triangular Toeplitz matrix generated by p/q 
satisfies the identity 

Tnip/q) =P V + R, P v £ e v , rank R < dcgp + 1 . 

Moreover, if the roots of q are real, then also the Hermitian Toeplitz matrix 
generated by $t(p/q) splits into 

T n QR(p/q)) = Q ip + R, Q V £G V , rank R < 2 rank R . 

Proof. By the fundamental theorem of algebra, any polynomial / : T — > C 
admits the splitting f(z) = ]X=i ( z ~ z i)- As a consequence, the rational 
function ^¥\ admits the simple fractions decomposition 

P(z) _ y^dogg Pi 



and Zi are the roots of q. 



where pi is the residual given by pi — (z — zi) ^\ 

Therefore, by the linearity of the Toeplitz operator T n : L 2 (T) — ^ M„ we have 
the identity 

Tn( P /q) = -Ei e A q %Zn(l/z k ). 

Now the existence of P v and R follows from Proposition 5.1. Finally, if all the 
roots of q are real, then the residuals pi are real, and 

T n (R(p/«)) = - EtY g- k (Kn(l/z k ) + I) , 
which together with Proposition 5.3 concludes the proof. □ 

Notice that, if all the roots of q are known, then one can explicitly compute 
the </?-circulant matrices of Theorem 5.4, obtaining not exactly the optimal but 
a rank bounded preconditioner for an important class of Toeplitz matrices. 

There follows one further lemma, whose proof can be found in [YR99] as 
underlined in [ZOT06]. 

For any e > there exist a, , bi such 
<ek~ a 
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Lemma 5.5 Let k £ {1, 


. . , n} and ael. Fo 


that 

p 




i=l 



with p < log e 1 (/3q + /?i log e l + @2 log n) , and the coefficients j3i depend only 
on a. 

By using the result stated in Lemma 5.1, we can reformulate the above 
Lemma 5.5. Let || • ||c denotes the Chebyshev norm on M,„ 



A — \Qij)i 



\M 



max |Oj 

ij 



Consider the lower triangular Toeplitz matrix T n = [{i — j) a ]i>j- Then, for 

(k) ~ 

any e > 0, there exist P^ ' e L and R^ of of rank one, such that 



\T n — 2^ik=l a kZn\t 



\C 



T 



V P 

Z^k = l ^1 



(k) 



Rk 



<e\\T n \ 



where p is bounded as in Lemma 5.5. By taking the transpose of T n and then 
subtracting by /, one immediately observes (via Proposition 5.3) that also the 
symmetric Toeplitz matrix T* — (\i — j\~ a )ij admits the decomposition 



\T~ — Cu, — R\\n < s\\T^ 



0(^ t \j(£ 



for a matrix R whose rank is bounded by 1p. 

From now on when referring to Lemma 5.5 we will always think at the latter 
two inequalities. 

Observe that the results obtained at this stage are enough to say that for 
any polynomial / and any symmetrized polynomial g(x) = /(|x|), the Toeplitz 
matrices T n (f) and T n (g) admit the decomposition 

T n = C v + R + E 

where C v € C v , R has sufficiently small rank and ||-E||c < £■ Notice furthermore 
that the same can be said for continuous symbol functions, since they can be 
approximated by polynomials (Weiestrass theorem). 

An even better result can be obtained for the Toeplitz matrix whose entries 
are positive integer powers of the indexes, namely 



LEMMA 5.6 Let T n = [(i- j) p ] t >j and T° = {\i-j\ p ) ir Then for any ip e T, we 
have T n = P^ + R and T° = Q v + R with P v , Q v e G v and rank fl < 2ranki? < 
2(p + 2). 



Proof. If we prove the decomposition for T n , the thesis for T£ follows because 
of the identity T* = T n + Tj. Fix p E N and consider the polynomial x v such 
that deg(Xip) < p + 1 and 



xA k ) -VX<p{k-n) = k p , 



,n — 1 



Call v p the first column of T n , i.e. T n 



«(X ¥ 



/ X V (0) \ 

x v (i) 

V X v ( n - !) J 



T n (v p , 0), then set 

Kx v . 



( X V (0) \ 

x v (-i) 



V X<p0--n) J 



2?, 



It follows that fJb(Xip) — n v a(xv) = ~n v v p , and therefore 

M<p(a(x<p),b{x>p)) =>K,p{v p ,Q), [T n (a(x<p) , Kx v )) , n >.p} = [T n (v p , 0), 77 v ] 
The thesis for T n follows by noting that rankT„(a(x v ), b(x<p)) < P + 2. 



□ 



Now consider a polynomial / of degree d By the previous Lemma we can 
affirm that the Toeplitz matrix whose entries are f(\i — j\) can be decomposed 
into the sum of a y-circulant matrix and a matrix R whose rank is bounded 
by ~^2 i= i i + 2 = 0(d 2 ). Such particular Toeplitz matrix is indeed a generalized 
KMS matrix. For the sake of completeness we recall that, given the matrix 

^(/,A) = (/(|i-i|)A^I 

a generalized KMS matrix is defined as 

m 
G n = £ ,JkK n (fk, Afc) 
fe=l 

where 7^ and Afc are all real and fk are polynomials of degree dk- 
It is clear that using both Proposition 5.3 and Lemma 5.5 we have 



\G n 



P« + R\\c < e\\G 



n\\C, 



*<£ ^ *^tn , 



with rank R = 0(log - (log - + log n) ^ 



k=i a k) 



Proposition 5.7 Set f(z) = log (2: - z ), z , zel Then, for any e > there 
exist P^ € C v and P £ with rankP e < loge" 1 (/3o + /3i loge^ 1 + /?2 logn) such 
that 

||T„(/)-P v -iy c < e ||T„(/)|| c . 

Proof. By the logarithmic singularity of / there follows the equality 

z k 
f(z)=1ogzo + J2u3;> 

k>l ° 



and thus for i > j, we have 

J-n\J)ij = 



(* - J>o 



(z-^r 1 ^^ 1 ) 



nj 



Note that the log zq term gives rice to a multiple of the identity, thus to an 
element of Q v . Therefore we do not care about it. By Lemma 5.1, the matrix 
Z n (zo) has the form P v + R for a P v € Q v and a rank one matrix R. Therefore 

(k) 

by Lemma 5.5, for any e > there exist P^ and a rank one Rk such that 



T n (f)-(P v + R)J2(P^ k) +Rk) 



< e\\T n (f)\\c ■ 



As a consequence we have the thesis, since (P v + R) 'J2 k (P^> + Rk) = P v + Re- 
with R F and P M as in the statement. □ 



24 



We can, finally, combine Theorem 5.4, Lemma 5.6 and Proposition 5.7, to 
obtain 

Theorem 5.8 Let / G £(T). For any e > there exist P v , Q v G G v such that 

||T„(/) -P v - R £ \\c < e\\T n {f)\\ c = 0(e) , 
\\T n ($f) -Q v - Re\\c < e\\T n ^Rf)\\ c = 0(e) 

with rank_R £ < 2rank_R £ < 2 log e -1 (a + feloge -1 + clogn) + d, and all the 
coefficients a, 6, c, d do not depend on n neither on e. 

The Hankel case can be discussed analogously. It is not difficult to check 
that H n (a,b) = JT n (Ja,b), for any a, b G C" with oi = b\. As before we call 
H(T n (a,b)) such matrix. Therefore we can reformulate the results obtained in 
this section simply multiplying them by J on the left, since this clearly does 
not affect the rank neither the arbitrariness of e. Moreover, we can write the 
Hankel matrix generated by the symbols 



C*(0) = i — — w' g( z ) = lo &( z - z o), z,z eT, 



l 

1 -~7*e ie 

in terms of H(Z n (X)). In fact, for instance, we have 

Hn(C^) = H n ~ l JZ n ( } r l ) = m"- 1 H(Z n (^ 1 )) (21) 

and a similar identity holds for H n (g) . 

We finally stress the fact that our initial problem is not well posed when 
/ is a rational function and the linear system is defined by the Hankel matrix 
with symbol /. In fact such matrix H n (f) has in general a small rank which 
equals the number of poles of / (due to the Kronecker theorem, 1881 [Kro81]) 
and therefore the linear system H n (f)x = y, when n is large enough, could even 
be unsolvable. 

Hartley-type algebras 

Let us conclude with few observations concerning Hartley- type algebras. In 
studying this case the arbitrariness of tp G T is crucial, in fact it allows us to 
use both circulant and skew-circulant matrices and thus to consider Hartley- 
type algebras. If "K is a generic Hartley-type algebra, recall that C* C "K with 

pe{-i,i}. ; 

Observe that by Proposition 5.3 we already know that the KMS matrix 
K n (\) admits the splitting K n (X) = H + R where ranki? = 2 and H is an 
element of C^ c "K. In fact we have K n (X) — C ip ((, ip (X)) + R where the matrix 
C V (S, V (X)), if G {1,-1}, is circulant or skew-circulant symmetric, respectively. 
Thus by the definitions in Table 4, when (p G {1, -1}, C v (£<p(\)) G "K. 

Let us summarize this remark into the following 

Lemma 5.9 Let A G K. For any Hartley-type algebra 'K there exists H G "K 
and R of rank two, such that K n (X) = H + R. 
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Proof. Specialize Proposition 5.3 for <p = 1 and <p = — 1 and use the definitions 
in Table 4. □ 

Nonetheless we would stress the fact that all the results we obtained in terms 
of symmetric </?-circulants may also be seen as involving Hartley- type algebras. 
One just need to specialize them for ip = I or ip = — 1. 

6 Conclusions 

We have tried to produce a first step towards the generalization of the ideas 
presented in [OT06] . After a brief overview about matrix algebras of low com- 
plexity, their generators and main properties, we have proposed a way to extend 
the applicability of the black-dot algorithm, proposed in [OT06] for the construc- 
tion of optimal rank circulant preconditioners for a Toeplitz system, to other 
types of linear systems (including Toeplitz plus Hankel like) and to precondi- 
tioners chosen in other low complexity matrix algebras. Then we have shown 
that, in fact, a suitable class of Toeplitz and Hankel matrices is indeed repre- 
sentable as the sum of a (p-circulant matrix and a small rank perturbation, for 
any tp of modulus one. Combining such representation for ip — 1 and ip = — 1 we 
then derive an analogous decomposition involving matrices from a Hartley-type 
algebra and a low rank perturbation. 

It is important to note that for a significant class of Toeplitz (and Hankel) 
matrices associated with a rational symbol, the optimal rank (^-circulant and 
Hartley-type preconditioners (as we called it) can be explicitly computed with- 
out the use of the black-dot method, provided that the symbol function and its 
poles are explicitly known. 
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